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1  INTRODUCTION 

This  report  discusses  the  work  accomplished  during  the  AFOSR  Grant 
AFOSR-84-0203,  entitled  "The  Influence  of  Electric  Current  on  Crack  Propagation  in 
Thermal  Fatigue  Tests."  The  goal  of  this  work  was  to  develop  a  method  for 
estimating  the  crack  tip  stresses  that  result  from  a  non-uniform  temperature 
distribution  near  the  crack  tip  that  occurs  when  resistance  heating  is  used  to 
thermally  cycle  fatigue  specimens. 

A  problem  associated  with  thermal  fatigue  tests  is  obtaining  sufficiently  rapid 
thermal  cycles  in  the  test  specimens.  In  order  to  overcome  this  problem,  large 
electric  currents  are  often  used  to  cyclically  heat  the  samples  since  this  mode  of 
heating  can  be  easily  controlled  and  rapidly  cycled.  However,  in  the  presence  of  an 
open  crack  an  otherwise  uniform  electric  field  is  disturbed,  causing  a  1  l^/r~ singularity 
in  the  current  which  results  in  localized  heating  at  the  crack  tip.  Thus,  cracked 
thermal  fatigue  specimens  which  are  electrically  heated  have  a  distinctly  non-uniform 
temperature  field  with  the  crack  tip  acting  as  a  heat  source  in  the  material. 

Recently,  experiments  were  conducted  through  a  program  sponsored  jointly  by  the 
Air  Force  Wright  Aeronautical  Research  Laboratories  and  the  Defense  Advanced 
Research  Projects  Agency.  The  aim  of  the  experiments  was  to  address  the  impact  of 
simultaneous  variations  in  temperature  and  stress  on  crack  growth  predictions  [1]. 
Temperature  variations  were  achieved  by  application  of  a  60  Hz  a.c.  current. 
Comparison  of  results  obtained  on  samples  which  were  thermomechanically  fatigued 
out-of-phase  (temperature  and  load  cycles  were  180°  out  of  phase)  to  those  fatigued 
in-phase  showed  a  significantly  greater  crack  growth  associated  with  the  in-phase 
testing,  while  no  discernible  difference  between  the  crack  growth  rates  of  out-of¬ 
phase  and  isothermally  tested  samples  was  observed.  An  additional  outcome  of  the 
testing  is  a  distinct  difference  in  crack  growth  behavior  between  specimens  tested 
out-of-phase  and  those  tested  either  isothermally  or  in-phase.  The  out-of-phase 
samples  exhibited  significant  shear  lips  while  those  tested  isothermally  and  in-phase 
exhibited  flat  crack  faces. 

A  possible  explanation  for  the  anomalous  crack  growth  behavior  is  that  thermal 
stresses  caused  by  the  resistance  heating  were  not  taken  into  account  in  calculations 
of  stress  intensities.  It  is  the  goal  of  this  work  to  estimate  the  contribution  of  the 
electric  field  to  the  stress  intensity  as  induced  by  the  associated  quasi-static 
temperature  field  to  estimate  the  magnitude  of  its  contribution  under  test  conditions. 
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The  problem  of  analytically  determining  the  mechanical  state  of  the  crack  tip  in  an 
electrically  heated  specimen  requires  the  solution  of  two  interrelated  problems: 
determining  the  disturbed  electric  field  in  the  crack  region  from  which  heat 
generation  can  be  calculated  and  secondly,  using  the  heat  generation  as  a  heat  source 
in  the  thermoelasticity  problem  to  calculate  the  stress  field. 


The  electric  field  problem  has  been  of  interest  in  fracture  mechanics  because  of  its 
role  in  the  widely  used  electric  potential  technique  for  monitoring  fatigue  crack 
growth.  Analytical,  experimental  and  numerical  procedures  have  been  developed  to 
determine  suitable  calibration  curves  which  relate  the  potential  drop  across  the  crack 
to  the  crack  length,  see,  for  example,  references  [2-6].  The  analytical  expressions 
generally  use  the  analogy  of  fluid  flow  around  a  flat  plate  and  a  mapping  technique 
such  as  the  Schwarz-Christoffel  technique  used  in  the  present  analysis  [7]. 


The  thermoelastic  problem  is  more  complex.  The  estimation  of  thermoelastic 
stresses  induced  by  a  heat  source  at  the  tip  of  a  crack  requires  determining  the 
solution  to  Poisson's  equation  to  obtain  the  stationary  temperature  field  and  the 
solution,  to  the  nonhomogeneous  biharmonic  equation  to  obtain  the  stresses  due  to 
the  temperature  field.  The  latter  part  uses  the  Duhamel-Neumann  analogy  and 
associated  modified  Airy  stress  function  to  separate  the  nonuniform  temperature 
distribution  from  the  elasticity  problem  [8-10]. 


The  first  work  dealing  with  the  two-dimensional  thermal  stress  fields  associated 
with  cracks  was  published  by  Florence  and  Goodier  [11]  and  Sih  [12].  The  method 
of  complex  representation  of  the  elastic  state  was  used  to  derive  the  theoretical 
stress  intensity  factors  in  infinite  two-dimensional  cracked  regions  for  both  the 
symmetric  and  antisymmetric  problem  classes.  Sih  gives  the  thermal  stress 
components  near  the  crack  tip  (for  the  geometry  in  Figure  1,  neglecting  higher  order 
terms)  as 
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where 


K  =  Re[  A] 

1  A 


(4) 


is  the  symmetric  stress  intensity  factor. 


-2 

K  =  =  lm[  A] 

2  A 


(5) 


is  the  antisymmetric  stress  intensity  factor  and  A  is  a  complex  constant  which 
depends  on  crack  configuration  and  the  form  of  the  prescribed  temperature  contraints. 


These  expressions  are  identical  to  the  general  elastic  crack-tip  stress  field  given  by 
Irwin  [13]  for  traction  boundary  condition  problems.  With  some  manipulation,  the 
corresponding  stress  distribution  in  terms  of  polar  coordinates  can  be  written  using 
identities  given  by  [14]  as: 
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where  K1  and  «2  are  as  defined  above. 


An  alternative  formulation  for  this  class  of  problems  was  presented  by  Bapu  Rao 
[15]  for  determining  the  state  of  thermal  streses  in  an  insulated  crack  in  an  infinite 
thin  plate  subjected  to  uniform  heat  flow.  In  this  case,  both  the  stationary 
temperature  and  stress  problems  were  formulated  in  the  elliptic  coordinate  system 
and  the  solution  developed  completely  in  elliptic  coordinates  using  biharmonic 
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functions  and  then  the  results  were  transformed  through  a  series  expansion  to  polar 
coordinates  with  the  crack  tip  as  the  origin.  The  final  results  obtained  for  the  local 
stresses  at  the  crack  tip  in  polar  coordinates  agree  completely  with  those  obtained 
by  Sih.  It  is  noted  that  the  real  parts  of  the  stress  field  are  identical  to  the  singular 
terms  of  the  Williams'  eigenfunction  relations  expressed  by  Torvik  [16,17]. 

Experiments  performed  by  Svoboda  [18]  showed  good  agreement  between  the 
stress  field  corresponding  to  uniform  heat  flow  presented  in  [12]  and  experimental 
data.  The  test  specimens  used  were  made  of  plexiglas  since  this  material  is 
essentially  brittle,  can  be  prepared  easily  and  is  relatively  isotropic.  In  this  case, 
the  complex  constant.  A,  corresponded  closely  to  that  obtained  by  Sih  for  the 
antisymmetric,  uniform  temperature  field  case: 

2 

iGa  a  sin  y 

A  =  - - - - —  •  VT  (9) 

(1  ♦  v  ) 
o 

where  a  is  the  thermal  expansion  coefficient,  v  is  an  elastic  constant  (v  =  3  -  4v, 

o  o 

plane  strain,  v  =  (3  -  v)/(  1  +  v),  plane  stress),  a  is  half  the  crack  length,  and  VT  is 
o 

the  undisturbed  temperature  gradient  directed  at  an  arbitrary  angle  y  with  respect  to 
the  crack  axis. 

Thus,  estimating  the  stress  field  induced  by  electric  current  flow  around  a  crack 
reduces  to  determining  the  coefficient  of  tie  singular  term  in  Williams’  eigenfunction 
series  caused  by  the  symmetric  temperature  field  induced  by  the  electric  current.  In 
addition,  the  solution  to  the  elasticity  problem  must  satisfy  appropriate  boundary 
conditions.  In  the  problem  considered  herein  the  boundaries  will  be  considered 
traction  free,  consequently,  any  stresses  that  are  calculated  exist  in  addition  to  those 
induced  by  mechanical  loading. 

Because  of  the  complicated  form  of  the  boundary  conditions  in  elasticity  problems, 
closed-form  solutions  to  problems  of  cracks  in  finite  regions  are  rarely  obtained. 
Several  investigators  have  studied  the  problem  of  thermoelastic  stresses  due  to 
various  temperature  fields  in  cracked  regions.  Konishi  and  Atsumi  [19]  used  the 
method  of  integral  equations  to  numerically  solve  the  antisymmetric  uniform  heat 
flow/insulated  line  crack  problem  for  a  semi-infinite  strip  with  the  crack  parallel  to 
the  surface.  Their  results  revealed  a  stress  intensity  factor  proportional  to  the 
temperature  gradient  times  the  crack  length  to  the  3/2  power,  consistent  with  Sih's 
result. 


Sekine  [20,21]  also  studied  the  two-dimensional  insulated  line  crack/uniform  heat 
flow  problem  in  a  semi-infinite  strip  with  the  crack  arbitrarily  oriented.  In  [20],  the 
temperature  was  prescribed  on  the  surface  of  the  crack  and  the  crack  was  modelled 
as  continuous  distributions  of  heat  sources  and  edge  dislocations.  Singular  integral 
equations  were  derived  and  the  solution  obtained  is  in  the  form  of  the  product  of 
the  series  of  Tchebycheff  polynomials  and  their  weight  functions.  In  [21],  the 
corresponding  problem  for  a  thermally  insulated  line  crack  was  solved  by  introducing 
a  continuous  distribution  of  sources  of  temperature  discontinuity  rather  than  that  of 
heat  sources.  Similarly,  Tchebycheff  polynomial  series  were  obtained.  The 

coefficients  of  the  polynomials  were  determined  by  solving  the  set  of  linear 
algebraic  equations.  The  stress  intensity  factors  for  the  symmetric  and 

antisymmetric  types  were  numerically  evaluated  and  compared  graphically.  The 
results  for  the  limiting  case  when  the  crack  is  situated  extremely  far  from  the 
surface  are.  again,  in  complete  agreement  with  Sih.  The  results  for  the  case  when 
the  uniform  heat  flow  is  perpendicular  to  the  plane  of  the  crack  are  in  agreement 
with  Konishi  and  Atsumi. 


Other  geometry  and  temperature  field  two-dimensional  crack  problems  have  been 
studied.  Sekine  [22]  also  studied  the  problem  of  a  heated  bounding  surface  and 
solved  it  numerically  again  using  singular  integral  equations  and  Tchebycheff 
polynomials.  Tweed  and  Lowe  [23]  used  transform  techniques  to  reduce  the  general 
problem  of  a  crack  in  a  half-plane  to  that  of  solving  integral  equations.  The  point 
source  problem,  was  solved  numerically  for  a  hollow  cylinder  of  infinite  length  by 
Herrmann  and  Kuemmerling  [24]  by  using  singular  integral  equations.  They 
emphasized  that  for  cracks  of  large  length  where  the  surfaces  strongly  influence  the 
stress  intensity,  this  numerical  method  is  inappropriate  and  suggested  using  finite- 
element  or  the  well  known  J-integral  (a  note  on  the  use  of  the  J-integral  for  thermal 
stress  problems  is  in  [25]). 


Finally.  Sumi  and  Katayama  [26]  studied  thermal  stresses  at  crack  tips  in  finite 
plates.  Their  results,  for  both  symmetric  and  antisymmetric  cases  were  obtained  by 
evaluating  the  unknown  coefficients  in  the  eigenfunction  series  using  a  modified 
mapping  collocation  method. 


The  role  of  electric  current  in  thermo-mechanical  fatigue  tests  has  previously  not 
been  established.  The  presence  of  the  crack  causes  a  singularity  in  the  electric 
current  which,  because  of  resistance  heating,  induces  an  even  more  severe  singularity 
in  the  heat  flux  at  the  crack  tip.  This  singular  heat  source  can  lead  to  significant 
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temperature  gradients  which  are  not  being  taken  into  consideration  by  temperature 
measurements  made  at  some  finite  distance  from  the  crack  tip.  Thus,  one  goal  of 
the  current  work  is  to  establish  the  magnitude  of  the  temperature  gradient  near  the 
crack  tip  in  order  to  estimate  errors  that  can  occur  in  measuring  crack  tip 
temperatures.  The  temperature  field  that  results  from  the  singular  heat  flux  also 
gives  rise  to  thermal  stresses  which  load  the  crack  tip.  These  stresses  will  affect 
how  the  crack  propagates.  Thus,  a  second  goal  of  this  work  is  to  establish  a 
method  for  estimating  the  crack  tip  stresses  that  are  induced  by  the  current/crack 
interaction. 

In  order  to  establish  these  estimates  of  crack  tip  temperature  and  stresses  in  a 
form  which  is  easily  applicable  to  a  variety  of  materials,  the  problem  is  simplified  in 
the  following  manner,  enabling  relatively  simple  closed  form  solutions.  The 
specimen  which  is  heated  electrically  may  have  the  edge  crack  geometry  of  Figure  2. 
The  problems  that  will  be  analyzed  will  be  concerned  with  the  circular  subregion  (R 

a 

of  Fig.  2)  of  radius  'a'  which  surrounds  the  crack  tip.  The  electric  potential  in  this 
region  is  assumed  to  be  the  same  as  that  which  occurs  for  a  crack  in  an  infinitely 
wide  plate.  The  current  terms  that  are  singular  are  extracted  from  the  infinite  plate 
solution  by  representing  the  solution  in  terms  of  the  eigenfunctions  associated  with 
the  cracked  circular  plate.  The  temperatures  in  the  cracked  circular  plate  are 
calculated  based  on  three  assumptions:  the  temperature  is  a  constant  reference  value 
on  the  outer  boundary,  the  crack  acts  as  a  perfect  insulator,  and  the  region  is  heated 
only  by  the  singular  heat  flux  term.  Once  the  temperature  field  is  determined,  the 
thermal  stress  fields  in  the  cracked  circular  plate  are  calculated  under  the  assumption 
that  its  boundaries  are  traction  free.  The  thermal  stress  problem  is  solved  using 
Williams'  stress  functions  which  satisfy  the  boundary  conditions  on  the  crack  face 
and  which  are  numerically  made  to  satisfy  the  boundary  conditions  on  the  outer 
boundary  in  a  least  squares  sense.  Only  the  lowest  order  stress  function  gives  rise 
to  singular  crack  tip  stresses  and  its  coefficient  is  proportional  to  K,  the  'stress 
intensity'  parameter  commonly  used  in  fracture  mechanics  to  assess  damage.  Thus, 
the  solution  results  in  a  relatively-  simple  expression  for  estimating  the  stress 
intensity  induced  by  local  heating  of  the  crack  tip  that  can  occur  while  using  electric 
currents  to  heat  thermal  fatigue  specimens. 

Stress  intensity  factors  for  some  nickel  alloys  are  calculated  using  the  data 
presented  in  [1],  and  are  compared  with  the  stress  intensity  factors  induced  by  an 
isothermal  mechanical  fatigue  specimen  of  the  same  type.  The  results  show  the 


thermal  stress  intensity  induced  by  the  electric  current  is  significant  and  is  additive 
and  therefore  can  contribute  to  an  increased  crack  growth  rate  in  electrically  heated 
thermal-mechanical  fatigue  tests. 

2  ELECTRIC  POTENTIAL  PROBLEM 

The  electric  potential  problem  is  solved  for  the  model  problem.  First  the  electrical 
potential  over  the  circular  region  is  established  through  a  mapping  technique.  From 
this,  the  boundary  condition,  E(a,0)  =  lid)  is  determined  for  the  circular  region 
containing  the  crack.  The  boundary  condition  is  then  used  to  numerically  evaluate 
the  Fourier  coefficients  in  an  eigenfunction  representation  of  the  electric  potential. 
This  approach  isolates  the  term  that  causes  a  singular  heat  flux  at  the  crack  tip.  It 
is  the  contribution  of  this  singular  heat  source  to  the  thermoelastic  stress  field  which 
will  be  evaluated  in  this  report. 

2.1  Formulation 

Two  regions  are  first  defined  (see  Figure  2).  R^  corresponds  to  the  circular  region 
surrounding  the  semi-crack  in  which  the  thermoelastic  problem  will  be  solved,  while 
Rq  is  the  complementary  region  of  the  plate,  not  including  R^, 

Ra  *  {(r,0)|O  <  r  <  a,  -*  <  9  <  n) 

W  L  L  / 

R  =  {(x,y)  1 0  <  x  y  <  -.  x,y  f  R  ) 

o  1  o  2  2  'a 

The  electric  potential  problem  is  then  stated:  determine  E(r,0)  such  that 

V2E(r,0)  =  0  on  R  (10) 

a 

and  satisfying  the  boundary  conditions 
F(a,0)  =  f  (9), 

where  f (9)  is  determined  below.  The  condition  of  perfect  insulation  (no  current  flow) 
across  the  crack  requires 


3E  3E 

_  -  0 


Finally,  the  potential  must  also  be  bounded  at  the  crack  tip. 
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|E(r,0)|  <  M  as  r  -)  0 


2.2  Determining  the  Boundary  Condition,  f(8) 

In  order  to  establish  the  appropriate  boundary  condition,  a  Schwarz-Christoffel 
conformal  mapping  is  used.  From  [7],  the  solution  to  the  problem  of  uniform  flow 
about  a  flat  plate  is  given  by  the  complex  potential 

H'(z)  =  \Mz2  +  a2)1/2  =  ^(x,y)  ♦  iy(x,y)  (13) 

here  ^(x,y)  is  the  velocity  potential,  fix, y)  is  the  stream  function,  Vq  is  the  far  field 

velocity  and  z  =  x  +  iy.  The  fluid  problem  is  analogous  to  the  electric  potential 

problem  where  V  is  equivalent  to  E  / L  (E  is  the  applied  voltage  and  L  the  specimen 
o  o  o 

length).  The  electric  potential  is  analogous  to  the  velocity  potential  which  is  the  real 
part  of  Viz).  Thus,  to  find  the  potential  on  the  boundary  r  =  a,  the  complex 
geometry  in  Figure  3  is  used: 


^(x,y)  =  Re[V]  =  ~^\PyP2  cos  '  ~ 


where  p ^  p  ,  6  ^  #2  are  defined  by  Figure  3.  On  the  boundary  r  =  a. 


f  I  =  —X1  5  +  4  cos#  cos  {  -  £  9  +  —  * 
,  L  '  2  2 


arc  cos 


-sin# 

5+4  cost 


=  m 


2.3  Solution 

The  solution  to  the  electric  potential  problem  is  expressed  in  terms  of  the 
eigenfunctions  of  the  problem  by  using  the  separation  of  variables  technique: 

E(r,#>  =  R(r)0(#>  (16) 


Problem  symmetry  leads  to 


E(r,#>  =  C  »Tc  rn/2  cos  n(#  -  n)l2 

n  n 


The  electric  potential  within  the  circular  region  is  obtained  by  determining  the 


.» f 
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coefficients  in  the  Fourier  cosine  series.  Since  the  functions  are  orthogonal,  the 
coefficients  can  be  evaluated  directly, 

a(2_n,/2  r  me  -  n) 

C  =  -  \  cos  -  f(0)d0  (18) 

"  n  J  2 

77 

by  performing  the  integration  numerically.  Results  of  this  numerical  integration  are 
tabulated  in  Table  1  for  values  of  n  up  to  32.  Neglecting  higher  order  terms,  the 
first  coefficient  (singular  current  coefficient)  gives  an  expression  for  potential  as  a 
function  of  ,r  and  6: 

-  E 

E(r,0)  2ar  -2  cos {6  -  n)l 2  (19) 

It  is  this  component  of  the  electric  potential  which  gives  rise  to  the  dominant 
(singular)  heat  source  at  the  crack  tip  and,  consequently,  it  is  the  contribution  of  this 
term  which  is  evaluated  in  this  report. 

3  THERMOELASTICITY  PROBLEM 

Once  the  potential  over  the  plate  is  obtained,  the  heat  flux  in  the  circular  region 
can  be  determined.  From  this,  the  thermoelasticity  problem,  which  includes  finding 
the  temperature  field,  is  formulated  using  the  Duhamel-Neumann  analogy  to  model  the 
problem  as  one  with  a  constant  temperature  field.  Using  this  approach,  an  Airy 
stress  function  formulation  may  be  used  in  which  the  temperature  field  produces 
surface  tractions  and  body  forces.  The  solution  to  the  resulting  non-homogeneous 
biharmonic  equation  is  determined  using  standard  techniques.  The  solution  to  the 
homogeneous  part  is  represented  in  terms  of  Williams'  eigenfunctions,  the 
coefficients  of  which  are  evaluated  numerically.  The  coefficient  of  the  lowest 
Williams'  eigenfunction  establishes  the  magnitude  of  the  singular  stress  field  at  the 
crack  tip,  and  it  is  this  quantity  which  we  seek  to  evaluate. 


3.1  Formulation 

A.  The  Temperature  Field 
The  heat  flux.  h(r,f?),  induced  by  the  current  is 

h(r,0)  =  <r(VE  •  VE) 


(20) 
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where  a  is  the  electrical  conductivity.  For  steady-state  temperature  fields,  the  local 
temperature  variation,  T(r,0),  is  related  to  the  heat  flux  by 

,  h(r,0) 

V2T(r,0)  - -  (21) 

k 

where  k  is  the  thermal  conductivity. 

The  temperature  field  must  satisfy  the  condition  of  a  constant  reference 
temperature  at  the  edge 

T  =  T  on  r  =  a,  -  n  <  6  <  n  (22) 

o 

perfect  insulation  across  the  crack  requires 


0T  3T 

_  .0.  o  <  r  <  a. 


finally,  the  temperature  must  remain  bounded  in  the  limit  as  r  approaches  zero, 


|T(r,0)|  <  M  as  r  -»  0. 


B.  The  Stress  Field 

Duhamel's  analogy  states  that  the  stress  field  in  a  body  subjected  to  a  temperature 
distribution,  T(r,0),  but  with  zero  body  force  and  zero  surface  tractions  is  the  same 
as  that  given  by  a  superposition  of  hydrostatic  pressure1 

oe--  pi.  rrd  -  0 

with  that  given  by  another  problem  in  which  the  same  body  is  kept  at  a  uniform 
reference  temperature  but  subjected  to  a  body  force  with  components 


0T 

F'  s  -  '  *  • 


To  satisfy  equilibrium,  the  hydrostatic  pressure  field  can  be  thought  of  as  having  body  force  components  F 
*  P  <0T/0r)  and  F^  *  ^(0t/0#)r  and  surface  tractions  having  a  normal  component  of  ~/#T.  Thus,  the 
superposition  with  the  second  problem  yields  the  original  conditions  (aero  body  force  and  surface  tractions) 
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V] 


m 


i 


m 

$ 


m 


mi 


m 

Mfcl 
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and  with  a  surface  traction  with  a  normal  component  of  /ST  where 


1  -  2v 


corresponds  to  plane  strain  and 


1  -  V 


to  plane  stress. 

The  solution  to  this  problem  is  solved  by  introducing  the  Airy  stress  function.  ♦, 
such  that 


1  84>  1  824> 

'  *  7  5T  *  ?  a?2  *  p 


a  i  a* 

r  <29) 

dr  r  dd 

Stresses  given  in  this  manner  automatically  satisfy  equilibrium.  Compatibility  is 
satisfied  if 

V4*t>  =  -  *V2T  (30) 


where 


1  -  V 


corresponds  to  plane  strain  and 


I'M 
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*  =  «E  (32) 

to  plane  stress. 

The  stresses  associated  with  the  Airy  stress  function  must  satisfy  the  modified 
traction  boundary  conditions: 

a r  =  /ST,  r ^  =  0  on  r  =  a,  -n  <  6  <  n  (33) 

o  q  -  /ST,  =  0  on  crack  faces,  6  -  ±  v,  0<r<a  (34) 


3.2  Solution 

A.  The  Temperature  Field 

Applying  Ohm's  law,  the  heat  flux  induced  by  the  current  is 

h(r,0)  -  a  —  (E  /L)2  (35) 

2r  o 

which  yields  the  axisymmetric  temperature  distribution  when  Poisson’s  equation  (21) 
is  solved: 

T(r,0)  =  -  -  (E  ID2  (a  -  r)  +  T  (36) 

k  2  °  ° 


B.  The  Stress  Field 

From  (30)  and  (36).  the  resulting  nonhomogeneous  biharmonic  equation  in  plane  strain 
becomes 


V4+ 


aEa 

2(1  -  v)kr 


(E  ID2 


(37) 


This  equation  is  solved  using  superpositioning  and  separation  of  variables.  The 
solution  to  this  equation  (for  details,  see  Appendix  A)  is 
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where  are  the  symmetric  Williams'  eigenfunctions  which  yield  the  relationships 


00  p  Q  0  ^ 

ZD  r<n_2)/2  L  (6-n)cos(n-2)  -  +  (n+2(-1)n)cos(n+2)  -  J 

n  2  2 

r=  1 


iVS 


00  p  ^  ^  "1 

ZD  r<n_2,/2  L  (n+2)cos(n-2)  -  -  (n+2(-1)")cos(n+2)  -  J 

n  2  2 


r=  1 


00  p  0  0 

r  d  =  y*  D  r<n_2)/2  L  (n-2)sin{n-2) - (n+2(-1)n)sin(n+2)  -  J 

r0  n  2  2 

r=  1  * 


(39) 


C.  The  Stress  Intensity  Factor,  K1 

Obtaining  a  stress  intensity  factor  requires  finding  the  coefficients  Dn  which  satisfy 
the  appropriate  boundary  conditions.  The  boundary  conditions  (see  Appendix  A) 
require 


00  r  9  9  ~\ 

ZD  a(n-2)/2  j.  (6-n)cos(n-2)  —  +  (n+2(-1)n)cos(n+2)  —  J 

n  2  2 

n=  1 


<rac  E  ,  aE 

-  (— )  - :  (3  ♦  2cos 9) 

2k  L  9(1  -  v) 


0°  p  ff  ff  -i 

y  Dn  a<n~2)/2  L  <n-2)sin(n-2)  -  -  (n+2(-1)n)sin(n+2)  -  J 

n=  1 


(40) 


<ra  E  ,  aE 

o\2  _ 


(— ) 


2k  L  9(1  -  v) 


(2  sin#) 


The  the  stress  intensity  factor,  K  ,  is  defined  by 


lim 


k,  -  ;^o\ 2- 


2 nr  a . j  evaluated  at  9  ~  0. 


(41) 
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From  (39). 

K1  =  4v27  Dr 
From  Appendix  A, 

a5/2  c  c 

-<ra  E  -  aE 

D,  =  — -  Mr  - - :  (3d,  -  2d,)  (42) 

1  2k  L  9(1  -  v)  1  2 

where  d1  and  d2  refer  to  the  normalized  coefficients  d^11  and  d^1*  which  correspond 
respectively  to  unit  pressure  and  sinusoidal  loading  on  the  boundary  r  *  a. 

Thus  the  stress  intensity  factor  can  be  written  generally  as 


*2  2  it  a  E<t  E  , 
—  -  (_0)2 

9  k(1  -  v)  L 


.5/2 


[3d, 


2d2] 


(43) 


for  the  plane  strain  case. 

The  generic  coefficients  d1  and  d2  were  determined  numerically  by  least  squares 
technique  (see  Appendix  B).  For  this  configuration, 

d  *  0.560764 
d2  *  -0.913788 


4  RESULTS 

The  key  results  of  the  benchmark  problem  are  the  following  two  simple  expressions 
which  give  estimates  of  the  error  in  the  steady  state  temperature  and  the  associated 
thermally  induced  stress  intensity  factor: 


T  -  T 


a  E  . 
—  (-2)“ 
2k  L 


aEa  E 

K.  *  0.231  -  (-2)' 

2k  L 


1 


(44) 


(45) 
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for  Poisson's  ratio  equal  to  0.3.  Applying  these  results  to  the  specific  data 
generated  by  the  thermal  mechanical  fatigue  tests  in  reference  1.  it  is  shown  that  the 
effect  of  the  current  on  the  local  temperature,  and  consequently  thermoelastic  stress 
intensity  factor,  can  be  quite  significant. 

The  materials  used  for  high-pressure  compressor  and  turbine  disks  and  spacers  in 
the  Air  Force  F100  engine  are  high  strength  nickel-base  alloys  IN100  and  Astroloy 
[1].  Fatigue  data  in  [1]  were  provided  for  IN100.  The  temperature  range  of  the 
tests  is  between  500  and  1200°F.  Materials  properties  data  for  IN100  are  listed  in 
Table  2  for  three  temperatures:  70,  600  and  1200°F.  Under  the  assumption  that  the 
sinusoidal  variations  of  heat  flux  (current)  and  temperature  are  in  phase  (steady  state 
thermal  analysis),  the  maximum  voltage  corresponds  to  the  maximum  temperature. 
Thus,  the  1200°F  data  were  used  for  calculating  the  magnitude  of  K. 

The  stress  intensity  factor  in  the  thermal-mechanical  fatigue  tests  at  termination, 
according  to  [1],  is  approximately  80  ksi  /"in?  A  significant  error,  then,  would  be 
about  10  ksi  /"in  (>  10%).  Using  the  above  expression  and  a  specimen  length  of  1.5 
in,  this  degree  of  error  would  require  an  applied  voltage  of  only  0.65  V  (for  a  critical 
crack  length  of  0.3  in.).  This  small  voltage  prompted  further  investigation  during 
which  it  was  learned  that  the  actual  voltage  drops  measured  across  the  test 
specimens  were  less  than  0.4  V  [27].  (Using  the  upper  limit  of  the  0  -  3  V  range 
quoted  in  [1]  would  yield  an  electric  current-induced  stress  intensity  in  this  material 
of  greater  than  150  ksi  /in".!)  Using  a  revised  voltage  range  of  0  -  0.4  V  and  the 
materials  data  in  Table  2,  the  most  severe  test  conditions  would  yield  a  current 
induced  stress  intensity  of  approximately  3.5  ksi  /in.  This  value  would  not  be 
expected  to  cause  significant  errors  in  the  experimental  methods  that  are  presently 
used  for  measuring  crack  growth. 

However,  the  expression  for  stress  intensity  factor  is  strongly  dependent  on  crack 
5/2 

length  (as  a  ).  Thus,  a  slight  extension  of  the  critical  crack  length  can  change  the 
previously  relatively  small  errors  in  K  to  significant  ones.  For  example,  extending 
the  critical  crack  length  from  0.3  in.  to  0.4  in.  results  in  an  approximate  10%  stress 
intensity  error  when  the  same  upper  limit  on  voltage  (0.36  V)  is  applied. 

The  simple  expression  for  current  induced  stress  intensity  presented  above  offers 
the  advantage  that  it  can  easily  be  applied  to  different  materials.  Due  to  variations 
in  materials  properties,  particularly  electrical  and  thermal  conductivity,  small  voltages 
may  have  pronounced  effects  on  the  electric  current  induced  K  for  other  materials. 
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Table  3  lists  several  commercially  available  high-strength  nickel  alloys  and  their 
compositions  for  comparison.  Some  of  these  alloys  are  used  in  gas  turbine 
manufacture.  Pertinent  materials  properties  for  each  of  these  alloys  are  listed  in 
Table  4.  The  room  temperature  and  600  °F  data  for  IN100  were  used  to  calculate 
stress  intensity  estimates  for  comparison  to  those  materials  listed  in  Table  4.  The 
ratios  of  K  to  the  appropriate  temperature  value  of  K  for  IN100  are  also  listed  in  the 
table.  It  is  noted  that  the  stress  intensity  estimates  of  all  materials  listed  are 
greater  than  that  of  IN100.  It  is  possible  that  the  resistivity  of  IN100  as  listed 
(calculated  from  data  obtained  during  the  fatigue  testing)  is  greater  than  the  actual 
value  due  to  equipment  configuration  and  testing  error.  This  would  suggest  that  the 
actual  stress  intensity  error  is  significantly  greater  than  the  5%  previously  stated, 
possibly  even  twice  as  high. 

Displacements  at  the  crack  edge  were  calculated  to  compare  with  test  observations. 
The  detailed  displacement  analysis  is  given  in  [35].  The  analysis  reveals  a  small 
crack  opening  displacement  of  approximately  0.0001  in.  These  results  are  consistent 
with  the  crack  opening  displacements  measured  on  test  specimens. 

5  DISCUSSION 

The  use  of  a  benchmark  or  model  problem  to  estimate  errors  in  stress  intensity 
factors  and  local  temperatures  in  thermal-mechanical  fatigue  testing  which  uses  large 
electric  currents  for  heating  has  shown  that  this  mode  of  testing  results  in  a  more 
severe  stress  state  than  if  conventional  heating  methods  are  used.  The  reason  for 
the  more  severe  stress  state  is  that  the  application  of  a  constant  voltage  across  a 
cracked  specimen  results  in  a  singular  current  and,  consequently,  a  singular  heat 
source.  The  singular  heat  source  does  not  result  in  a  temperature  singularity, 
however,  so  the  temperature  remains  bounded  at  the  crack  tip. 

The  relative  discrepancy  between  conventional  heating  fatigue  testing  and  resistance 
heating  fatigue  is  quite  dependent  on  crack  size  (as  a  )  and  materials  properties: 
primarily  thermal  and  electrical  conductivity  and  linear  thermal  expansion.  Materials 
exhibiting  low  electrical  resistance  and  high  thermal  expansion  are,  in  general,  more 
susceptible  to  this  mode  of  thermal  loading. 

Simple  expressions  for  estimating  the  contribution  of  electric  current  to  the  stress 
intensity  factor  and  local  temperature  field  were  developed.  These  expressions  are 
estimates  in  light  of  some  basic  assumptions  made  during  problem  formulation.  The 
validity  and  consequence  of  some  of  these  assumptions  will  be  discussed  here. 


-f.  -  .  -C.  ■r.  r-.  * .  J~.  ■r.  -',-"7!  V  IT,  -T.  ,  ■r.  --  -- 
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Beginning  with  the  electric  potential  problem,  we  first  assumed  that  skin  effects 
(the  concentration  of  electric  current  near  surfaces)  can  be  neglected.  This  allows 
the  problem  to  be  configurated  as  one  in  two-dimensional  plane  strain.  In  reality, 
however,  the  skin  effect  is  known  to  be  important  in  some  cases.  If  the  skin  effect 
is  important  in  the  present  system  the  consequence  is  a  hotter  surface  due  to  the 
higher  current  density  in  that  area.  One  example  where  the  skin  effect  is  important 
is  in  a  cracked  plate  with  uniform  applied  magnetic  field  [34].  The  increase  in 
current  density  uniformly  from  the  center  to  the  outer  edge  in  this  case  was 
determined  by  finite  element  to  be  approximately  20%.  It  is  possible  that  this  edge 
effect  may  be  related  to  the  shear  lip  formation  observed  in  some  of  the  fatigue 
specimens  in  that  K  could  be  higher  near  the  surfaces  of  the  specimens. 

Secondly,  the  assumption  that  higher  order  terms  in  the  electric  potential  series 
expression  can  be  neglected  when  calculating  the  heat  source  induces  an  error  when 
calculating  the  stress  intensity  factor.  The  trouble  with  including  more  than  a  limited 
number  of  terms  in  the  series  is  that  the  problem  would  then  have  to  be  solved 
completely  numerically.  By  using  the  current  approach,  relatively  simple  expressions 
were  derived  for  the  temperature  error  and  the  current  induced  stress  intensity  which 
may  be  easily  evaluated  for  any  material.  However,  these  results  are  not  exact  and 
only  provide  an  estimate  of  the  magnitude  of  this  effect. 

It  is  noted  again  that  although  the  greater  crack  growth  rates  are  observed  in  the 
in-phase  testing,  it  is  the  out-of-phase  test  samples  that  exhibit  the  anomalous  crack 
growth  behavior.  Time  histories  of  temperature  and  electric  current  that  were  taken 
during  a  crack  growth  test  [27]  are  plotted  in  Figure  5.  It  is  apparent  from  this  data 
that  the  assumption  that  the  current  and  the  temperature  are  in  phase  is  inappropriate 
under  actual  test  conditions.  In  fact,  it  is  clear  from  the  phase  lag  between  current 
and  temperature  that  thermal  transients  play  an  important  role  in  these  cases. 
Moreover,  the  temperature  information  in  the  figure  may  be  used  to  estimate  the 
relative  magnitude  of  the  time  dependent  term  in  the  heat  equation.  Averaging  over 
a  circle  of  radius  a,  the  heat  flux  generated  by  the  electric  current  (15  °F/sec)  is  of 
comparable  magnitude  to  the  thermal  transient  term  (12  °F/sec).  This  indicates  that 
the  assumption  of  steady  state  in  the  thermal  analysis  should  be  reconsidered  and 
that  transient  thermal  effects  should  be  taken  into  consideration  in  order  to 
accurately  relate  analyses  to  the  laboratory  situation. 

In  conclusion,  a  first  estimate  to  the  role  of  electric  current  in  thermal  fatigue  tests 
has  been  obtained.  It  was  calculated  under  the  assumption  that  the  singular  heat  flux 


term  would  dominate  the  near  field  behavior  at  the  crack  tip.  Appropriate  linear 
differential  equations  were  solved  analytically  up  to  the  point  in  the  thermoeiastic 
problem  where  the  boundary  conditions  had  to  be  satisfied.  A  numerical  technique 
based  on  a  least  squares  approach  was  used  to  satisfy  these  constraints. 


Since  little  work  has  been  done  in  this  area  the  resulting  solution  is  fundamental  in 
that  it  may  be  easily  used  to  establish  thermal  stresses  for  a  variety  of  test 
configurations.  Consequently,  it  may  be  readily  utilized  to  assess  the  validity  of 
thermal  mechanical  testing  procedures  which  incorporate  resistance  heating  and 
neglect  thermoeiastic  stresses. 
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Table  1 

Electric  Potential  Coefficient  Values 


Coefficient  Number 
(n) 


Value 


.970527300 

1.414210000 

.341726000 

-.353554000 

-.341728200 

-.044193700 

.069035320 

-.011048700 

-.054572040 

-.003452562 

.028407400 

-.001208524 

-.021641700 

-.000453151 

.014980200 

-.000178091 

-.011750900 

-.000072299 

.009123044 

-.000030232 

-.007432150 

-.000012750 

.006108524 

-.000005576 

-.005137681 

-.000002395 

.004368960 

-.000001079 

-.003766370 

-.000000430 

.003277960 

-.000000251 

-.002879890 


Notes: 


I', 


n [6  -  v ) 


C  =  A  —  a(2"n)/2  for  n  =  0.1, 2, 3,... 

n  n  L 


E  a  E 


C  =  0.9705273  -2-,  -2  in  volts/in;  C,  =J~2a  -2 
o  L  L  ’  V  l 


*  -  •  -  j*  -  -  -  •  W  r  i  .  -  h.  -  ->  -  -  -  -  -■  -  -*  -  •  -  •  -  -  w  *  n  .  *  /v  -  »*L\  -N  _  AN  -N  _  *  lN 
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Table  2 

Properties  of  IN  100 


Temperature 

<°F) 

Linear  Thermal  Expansion 
Coefficient,  from  70° 

(10-6  in/in/°F) 

Thermal 
Conductivity 
Btu-in/hr-ft  -  F) 

Electrical 
Resistivity 
(fiQ- cm) 

Young's 
Modulus 
(106  psi) 

RT 

6.811'21  ' 

80O) 

195(2'3) 

29.5* 11 

7.5(1) 

100(1) 

208(3! 

27.0(3) 

1200 

8.0(1) 

130!1) 

220l2,3) 

24.0(2-3) 

Notes: 


m 

1.  Aeronautical  Vest-pocket 

*  o 

p.  106-107. 

*■ 

.«  V 

2.  Extrapolated. 

Handbook,  P&WA  Publication.  July  1979, 


» 

*: 


3.  Calculated  from  data  obtained  from  J.  Warren. 


Table  3 

Major  Constituents  of  High  Strength  Nickel  Alloys 


Material 


Other 


Reference 


IN100 
Nickel 
K  Monel 
Monel 
Hastelloy  B 
Hastelloy  C 
"Inconel" 
Inconel  600 
Inconel  718 
Inconel  X750 
Incoloy  800 


55 

100 

63-70 

63-70 

60 

60 

72 

(Ni+Co>=72 

55 

(Ni+Co)=70 

33 


(3Mo,4  Ti,5AI) 


(~  25  Cu) 
30  Cu) 
(28  Mo) 
(18Mo,4W) 


(5  Nb) 


1 

28 

28.29 

28.29 
28-30 
28-31 
28-31 
32 
32 
32 
32 


r 

i1 


Table  4 

Materials  Properties  of  Several  Nickel-Base  Alloys 


I 


Linear 

Thermal 


Material 

Expansion 
Coefficient 
(10-6  in/in/°F) 

Electrical 
Resistivity 
(10-6  O-cm) 

Thermal 
Conductivity 
(Btu-in/in-ft  -  F) 

Young's 
Modulus 
(106  psi) 

Temp 

(°F) 

Ni 

7.4 

6.84 

637 

30 

70 

K  Monel 

7.8 

58.3 

130 

26 

70 

Monel 

7.8 

48.2 

173 

26 

70 

Hastelloy  B 

5.5 

135 

78.2 

28.5 

70 

Hastelloy  C 

6.3 

133 

86.9 

28.5 

70 

Inconel 

6.4 

98.1 

104 

31 

70 

Inconel  600 

5.8 

103 

103 

30 

70 

7.9 

107 

133 

27.5 

600 

8.6 

113 

172 

24.8 

1200 

Inconel  718 

7.8 

125 

77 

29 

70 

8.0 

129 

111 

26.7 

600 

8.6 

134 

147 

23.7 

1200 

Inconel  X750 

6.9 

121 

83 

31 

70 

7.5 

126 

109 

(28.5) 

600 

8.4 

131 

143 

25.5 

1200 

Incoloy  800 

8.7 

98.8 

80 

28.4 

70 

9.0 

111 

115 

25.4 

600 

9.6 

124 

152 

22.3 

1200 

K  refers  to  the  stress  intensity  factor  calculated  for  INI 00. 
IN 


*■'  to 


*  »'*'  .*» 
s-C  ju  J>V  j 


**-  ^.y. 


K/K. 


3.9 

2.0 

1.89 

1.16 

1.21 

1.51 

1.28 

1.57 

1.6 

1.83 

1.53 

1.54 
1.66 
1.6 
1.2 
2.4 
1.8 
1.6 


'^  r  , 
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Figure  1.  Thermal  Stress  Component  Geometry  Near  Crack 


a.  Mapping  Geometry 


b.  Problem  Geometry 


ft 


Figure  3.  Relation  of  Schwarz-Christoffel  Mapping  Geometry  to 
the  Problem  Geometry 


Heat  as  a  Function  of  Position  for  Various  Numbers 
of  Terms  in  the  Fourier  Cosine  Series 


\ 


Temperature  (°F) 


\ 


Current  (amps) 
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Appendix  A 

DETAILS  OF  THE  THERMOELASTIC  SOLUTION 


From 


A  ~aE  - 

v4+  =  - - :  V2T 


<1  -  v) 


and 


V2T  = 


a  E  _  a 

o\2 


- (-2)“  - 

2k  'l  r  ' 


the  governing  equation  can  be  written  as 


V4*  = 


a  a  E 


(E  ID 


2  a  6 


2k(1  -  v)  o  r  r 


where  ♦  must  satisfy  the  traction  conditions  of  (33)  and  (34).  Using  the  principle  of 
linear  superposition,  let 


<t>  =  4>  +  <t> 

H  P 


where  4>p  is  the  particular  solution  satisfying 


B 


and  ♦  is  the  homogeneous  solution  satisfying 

H 


v  *h  "  °- 


A  solution  for  4>p  is 


* 


Br' 


P  9 


The  homogeneous  solution,  ♦  can  also  be  broken  into  component  parts: 

H 


.  wvc-i 

m 
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M 


4>  =  <i>  -f  <t> 

H  HI  H2 

where  4>  is  the  Williams'  eigenfunction  solution  (♦  automatically  satisfies  <34)  and 
its  coefficients  are  numerically  determined  such  that  (33)  is  satisfied  in  a  least 
squares  sense)  while  frees  up  the  crack  faces  from  the  tractions  imposed  by  ♦  . 

Hi  p 

By  separation  of  variables,  an  expression  for  <t>  is: 

H  1 


Br 

*H1  =  T  cos* 


It  is  noted  that  the  combination  ♦  +  satisfies  the  crack  face  conditions: 

p  HI 


a  q  =  *  yST  =  /3T  at  8  =  ±  n 


3  ,  1  3$  x 

—  )=0  at  0  =  ± 
dr  \  r  30  / 


The  appropriate  boundary  conditions  for  <t>H2  are  derived  from 


|  1  d  id2 

al  'UVa?11*-*  *h.  *  '/rr-zn 


or.  after  incorporating  the  above  expressions  for  ♦  and  ♦  and  evaluating  at  r  =  a, 

P  H 


I  !5 H2  +  1  5  *H2  I 
r  3r  7  dff 7  '=» 


Ba 

—  (3  +  2cos0) 
9 


b)  r  q  I 
rfr  r=a 


a  i  a 

- r - (<f>  +  <t>  +  <t>  )]  =  0 

3r  r  d6  p  Hi  H2 


a  i  a*  ,  2 

■  f,  ‘r  V  ■  •  i Ba  sin* 


mm. 


V. 

I'J 

«y 

« 

w 

‘XI 

1% 


y 


Thus,  the  Williams'  eigenfunction  solution  must  satisfy  boundary  conditions  of  the 
form: 


a  1  = - (3  +  2cos#) 

r  r=a  g 


Ba 

= - 2  sin# 


Tr#  r=a  g 


In  order  to  generalize  and  further  simplify  the  numerical  evaluation  of  coefficients, 
D  ,  in  the  eigenfunction  solutions,  ♦LJ_,  these  coefficients  can  be  expressed  as  the 
sum  of  normalized  coefficients  d^  and  d2<n)  which  correspond  to  unit  pressure  and 
sinusoidal  loadings,  respectively,  on  the  boundary  r  =  a.  In  other  words,  if  on  the 
boundary  r  =  a, 

00  00 

K>  ■  Z  diW  ■  ’•  <',»>,  *  Z  d,"”  r„'«  ■  »• 

1  n-1  1  n=1 

and 

00  00 

(<r  )  *  y*  d  <n>  X  (#)  =  cos#,  (r  J  =  y*  d,(n)  V  (#)  =  sin#, 

r  _  2  n  r  O  _  /  ■  2  7  n 


2 


2  n=  1 


For  XJ8)  and  yn<#)  representing  the  Williams'  eigenfunction  series  for  and  r^, 
respectively,  the  coefficients  Dn  in  (40)  are  such  that 


=  -  —  (3d  ♦  2d  n  =  1,2 . oo. 

9  1  2 


Numerical  calculation  of  d^"1  and  d2<n)  is  discussed  in  Appendix  B. 

— 1/2 

Substituting  for  B  and  dividing  by  a  ,  the  singular  coefficient,  becomes 


D,  -  -  -  »5'2  Ml2  (3d,"1  •  2d 

1  2k  L  9(1  -  v)  1  2 


/■  >V  j- 'i,j V-  3-  r\'" 


mnuiviui'j'u 


U«'J  »'J»U 


r  i  *"j  vj  r:  v*  ■ 
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Appendix  B 

METHOD  OF  LEAST  SQUARES:  CALCULATIONS  AND  CONVERGENCE  CHECKS 


The  method  of  least  squares  was  used  to  determine  the  coefficients  d^"1  and  d2<n> 
so  that  the  appropriate  boundary  conditions  presented  in  Appendix  A  are  satisfied. 
The  least  squares  method  constructs  a  series  representation  which  fits  the  boundary 
conditions  over  the  interval  of  interest  with  a  minimum  of  error.  If  the  boundary 
conditions  are  written 

a  -  f(0),  r  . .  =  g{0)  on  r  =  a, 

r  r(J 

3 

the  error  involved  in  summing  a  ,  and  r^,  at  the  boundary  r  =  a  is 

it 

n  00  00 

<  ■  J  {[  «»>  -  Z  d-r  “»p  *  [  -  Z  V  «>?}  as 

0  n= 1  n=  1 

where  f ($)  and  yJ8)  represent  the  Williams’  eigenfunction  series  for  and  r  q, 
respectively. 

In  order  to  minimize  t  with  respect  to  each  coefficient  d^,  take 

00 

'  X  V«W)]  rJe) 

n=  1 

1  y  (0)  }  6d 

A  ? 

=  /  d  \  (y  ^  *  y  y 

/  -  .  n  J  7  nr  m  '  n '  m 
n=  1  0 


£-ft-t 


m 


00 

-  2  [  -  Z  dn> 


n=  1 


or. 


n 

s 


[f (d)f  *  g(d)y  (eme 


If  the  array  I  is  defined  as 


The  calculations  must  be  performed  separately  for  the  unit  pressure  and  sinusoidal  loading  conditions.  The 
method  in  the  text  describes  the  general  p^cjcedure  to  be  followed  in  both  cases.  For  ease  of  notation,  we 
have  used  d^  in  place  of  either  d  ^  or  d  in  this  appendix. 


a"* 


.  -\r«  T,  /.  -rm  .  ,  /.  /.  v.  -  .-  Tip:  v,  v,  s  wv 

*-“*-*■  "jjjAiAa  '-a  jjuJi  iJaAk  flaSiaa  fU  nfcha  rJakafri*  ar.afU  t 


and  the  constraint  vector  components  as 


W 


(f (9)f  +  g(0)y  )d0 


Then  the  error  is  minimized  in  a  least  squares  sense  if 

N 

Zd  I  =  F  as  N  4  oo 

n  nm  m 

n=  1 

This  equation  can  be  solved  for  the  unknown  coefficients  dn  for  various  values  of  N 
using  a  standard  LU-decomposition  method. 


To  check  the  numerical  convergence,  three  tests  were  performed: 

1.  Stresses  and  errors  were  evaluated  when  functions  f  *  1,  g  *  sin20  were 

inputted.  Results  for  the  cases  N  =  4  and  N  =  16  are  plotted  in  Figures 
B-1  and  B-2  along  with  the  known  functions  for  comparison.  Excellent 
correspondence  between  the  known  functions  and  the  series 

representations  for  N  =  16  is  noted,  while  the  case  N  =  4  still  results  in 
rather  large  error  values. 

2.  The  coefficients  were  calculated  for  a  test  case  in  which  the  exact 
coefficients  were  known.  These  were  used  to  generate  boundary 
conditions  to  see  if  the  same  coefficients  were  rendered.  Results  are 

-5 

tabulated  in  Table  B-1  and  show  error  values  on  the  order  of  10  . 

3.  The  coefficients  were  used  to  compute  the  boundary  stresses  to  compare 

with  the  theoretical  boundary  conditions.  Figures  B-3  and  B-4  show  the 
results  for  four  values  of  N.  The  results  show  that  the  boundary 
conditions  are  satisfied  quite  well  for  N  as  small  as  16.  Excellent 

agreement  is  noted  for  higher  values  of  N. 


The  convergence  checks  illustrate  quite  rapid  convergence  and  excellent  behavior  of 
the  functions  with  approximately  zero  error  accomplished  with  as  low  as  32 
coefficients. 


Figure  B-3.  Solution  Check:  o  at  Boundary  r  =  a  for  Various  Values  of  N 


